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Abstract —We present a longitudinal MR simulation framework 
which simulates the future neurodegenerative progression by 
outputting the predicted follow-up MR image and the voxel- 
based morphometry (VBM) map. This framework expects the 
patients to have at least 2 historical MR images available. The 
longitudinal and cross-sectional VBM maps are extracted to 
measure the affinity between the target subject and the template 
subjects collected for simulation. Then the follow-up simulation 
is performed by resampling the latest available target MR image 
with a weighted sum of non-linear transformations derived from 
the best-matched templates. The leave-one-out strategy was used 
to compare different simulation methods. Compared to the state- 
of-the-art voxel-based method, our proposed morphometry-based 
simulation achieves better accuracy in most cases. 
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I. Introduction 

Neurodegenerative diseases are a category of progressive 
loss of structure or function of neurons, including Alzheimer’s 
Disease (AD), Parkinson’s disease (PD), Huntington’s disease 
(HD) and Spinal muscular atrophy (SMA), etc. Neuroimaging 
has advanced the analysis of neurodegenerative processes pro¬ 
foundly in the past two decades with large-scale group studies 
E), IT E) The large-scale neuroimaging computing methods 
of neurodegenerative diseases can be generally categorised into 
cross-sectional and longitudinal. Though the majority of the 
recent neuroimaging studies focused on the cross-sectional 
group comparison with regional measurements m, [51, m, 
□ the longitudinal analysis of brain tissue changes is effective 
in describing an evolving neurodegenerative process ©, 0, 
EB. A number of investigations have attempted to analyse 
the longitudinal changes from neuroimaging biomarkers at 
serial time-points and to make predictions of the underlying 
process ED, ED, ED- However, it is difficult to evaluate 
any such predictions in practice due to the lack of ground 
truth. We suggest that simulation of follow-up brain images 
can be useful for validating the predictions obtained from 
neuroimaging biomarkers. 

The longitudinal simulation of the follow-up structural mag¬ 
netic resonance (sMR) images can be performed by resam¬ 
pling the baseline sMR image with a weighted average of 
the longitudinal non-linear transformations from a template 
population 031 • When the template population is sufficiently 
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large to include a majority of neurodegenerative changes, the 
confidence in simulated results would heavily depend on the 
metrics used to measure the affinity between the target patient 
and the template population. 

Previous frameworks of longitudinal sMR simulation gen¬ 
erally depended on two assumptions about the progress of 
neurodegeneration: (1) a common rate of atrophy is shared 
across all subjects with the same diagnostic label; (2) brains 
with similar morphology evolve in a similar way ED, LET 
m. Alzheimer’s disease (AD), which is the most common 
neurodegenerative disorder, progresses in a special pattern 
which has multiple stages. The neurofibrillary tangles (NFT), 
which are thought to contribute to local atrophy, spread from 
memory related areas towards areas in the medial temporal 
lobe, the parietal cortex and the prefrontal cortex fl6l . Thus, 
patients in different stages of progression stages might not 
have similar atrophy rates across all brain regions. There 
might also be difficulty in comparing the local structural 
morphology accurately between different patients only based 
on the sMR intensities due to the inter-subject variance in the 
original structural appearances. Besides the sMR intensities, 
the historical progression of the same patient (longitudinal) as 
well as the difference between the current state of the patient 
and a normal brain template may generate the simulation from 
an alternative perspective. 

In this study, we present a proof-of-concept framework to 
simulate neurodegeneration in follow-up with sMR data. An 
overview of the framework is illustrated in Fig. [T] We hy¬ 
pothesise that brains with similar detected cross-sectional and 
longitudinal atrophic deformations would have similar follow¬ 
up evolution. The cross-sectional changes are derived by a 
symmetric non-linear registration from a standard template to 
the subjects; the longitudinal registration between two serial 
sMR images of the same patient is used to extract the intra¬ 
subject serial changes. We apply the voxel-based morphometry 
(VBM) to measure both types of structural brain deformations. 
Each subject recruited in the template population is expected 
to have at least 3 serial sMR visits. Both the longitudinal 
VBM map and the cross-sectional VBM map are collected 
for each template subject. The simulation of the future sMR 
volume of a target patient is performed by resampling the 
latest available sMR volume with an average weighted sum 
of the longitudinal transformations in the following period in 
the template population. Only a minority of the best-matched 
templates in both VBM maps are selected to contribute to the 
simulation results. The cut-off threshold depends on the size 
of the recruited template population. 
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Fig. 1. An illustration of the proposed framework of the longitudinal MR simulation. 


II. Methods 

A. Preprocessing 

Each template subject Si is expected to have at least 3 
serial sMR visits taken at time-points a, 6, c and a < b < 
c A b — a ~ c — b. Each template MR volume is skull- 
stripped using the Brain Extraction Tool with a standard 
space pre-masking applied 03 . Then all template MR images 
are affine registered to the MNI152 template brain space 
with FLIRT lfT8l . A symmetric diffeometric registration is 
performed between each pair of the standardised MR images 
[. I - a ^, 1-^] and [if ^, I- c ^] to obtain the longitudinal brain tissue 
deformations X j ab ^ and T^ bc \ The symmetric registration 
firstly registers the two adjacent MR volumes to a mid-way 
space to ensure the inverse consistency of the forward and 
backward transformations which was known to be important 
for atrophy calculation [19), [20]]. The MNI152 template M 
is symmetrically registered to each l\ b ^ image to obtain the 
cross-sectional deformation r p^ Mb \ 


B. Template Weighting with Voxel-Based Morphometry 


The voxel-based morphometry (VBM) maps [2T), l22l were 
calculated on X j ab ^ and Xrespectively as the Jacobian 
determinant of the spatial transforms. U is the displacement 
of X. The displacement tensor of U over time t is represented 
as 
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The determinant det(^(x,t)) = det(I + VC/) is calcu¬ 
lated at each voxel and forms the VBMs. To measure the 


distance between two detected tissue deformation Ti and Tj, 
the squared Euclidean distance of the VBMs is calculated 
within the area of a dilated brain mask p in the standardised 
space as D(Ji,Jj) = - di j ^) 2 /\p\, where d^ is 

the Jacobian determinant or Ti at voxel v. The area near the 
outer boundary of the dilated brain mask is used to capture 
the deformation on or near the cortical surface. 


C. Follow-Up MR Simulation 

To simulate the follow-up MR volume of the target 
patient S*, we require the patient to have two MR volumes 
available and 7* and the interval (b — a) to approximately 
equal the intervals in the templates. The symmetric longitudi¬ 
nal and cross-sectional deformations are obtained as T£ b and 
T^ b . The VBMs of both transformations are correspondingly 
computed as J% b and J^f b which are used for collecting the 
masked longitudinal distances X>- on ^(*, i) = D(J% b , Jf 6 ) and 
the cross-sectional distances D^ ross {~k,i) = D{J^ b , X Z M6 ). 
Two types of collected distances are respectively normalised 
within [-1, 1] as D = (D — D) / max(D — D) where D is the 
mean value of D. Then both distances are summed with rela¬ 
tive weights to obtain the combined D[ = aD l ° n9 + /3f)9 ross ? 
a + (3 = 1. k nearest neighbours of S'* are selected to form 
a new template set K. The majority of the templates are 
dropped out for one simulation because the distant subjects 
would introduce bias rather than contributing to the simulation 
accuracy, even Gaussian distributed weights are used. Based on 
D[, the follow-up deformation T b c of subject S* is computed 
as an average weighted sum of the template transformations 
T bc as 

^... E, £W ir X e- D!/ » 
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TABLE I. The visual check of a successful simulation and a failed simulation. The column Year 1 is the MR image at the second 
time-point; the column Year 2 is the real follow-up MR image of this patient and the column Simulated is the predicted Year 2 MR 
image. The YBM maps are overlaid on the Year 2 and Simulated images. The failure was probably caused by the spatial mismatch 

INTRODUCED IN THE AUTOMATIC AFFINE NORMALISATION. 



where T?* is the forward transformation from the template 
image If to the target image If; g is the Gaussian kernel 
density which is set to 0.5 m . The follow-up image If can be 
simulated by resampling I b with Tfc. Notably, in this proposed 
framework, subjects with different diagnostic labels are not 
required to be computed separately, since there might not 
be clear boundaries between the atrophic progression patterns 
from different diagnostic labels. For example, some follow¬ 
up transformations of early AD patients may contribute to 
simulate the future neurodegeneration of a late MCI patient. 
With the drop-out threshold \K\, the distant subjects are 
expected to be filtered out before the template merging. 

III. Evaluation and Results 

We recruited Grad Warped MR images with N3 Correction 
from the publicly available ADNI 1 and ADNI GO dataset 
(http://adni.loni.usc.edu/) [23 1. The slice thickness of all MR 
volumes is 1.3mm. We kept 60 subjects who had at least three 


continuous MR visits with an interval of 12 months available, 
resulting in 180 MR volumes and 120 MR longitudinal pairs 
in total to be registered. After transforming each MR volume 
to the standard MNI152 space, skull stripping was conducted 
with a fraction threshold of 0.3. We used the MNI152 template 
with 2mm slice thickness for the image normalisation. 

We used the leave-one-out strategy to evaluate the frame¬ 
work. The I a and I b of each selected patient were used as 
the inputs to simulate the unknown follow-up I c . The rest of 
the patients were used as the simulation templates. All the 
transformations and images belong to the testing target are 
excluded from the template construction. The simulated I c is 
then compared with the real follow-up image as well as the 
image derived by registering I b to the real follow-up image. 

The symmetric non-linear registration was implemented 
with the Advanced Normalization Tools (ANTS) (24) . The 
integration of the entire framework was implemented based 
on the Nipype framework (25). 
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Fig. 2. The plot of the three types of MR image distances considered for this evaluation: The distances between the predicted images and the real follow-ups 
(P-B); The distances between the predicted images and the registered 1 year MR images (P-rA); The distances between the registered 1 year MR images and 
the real follow-up MR images (REAL) 


Examples of a successful simulation on an AD patient and a 
failed simulation of a normal control patient are presented in 
Table. [I] In each case, the Year 1 MR image ( I b ) is shown 
as a reference; The real follow-up Year 2 MR image (/ c ) 
and the simulated follow-up MR image (/*) are overlaid on 
the VBM map extracted between them and the Year 1 MR 
image to visualise the longitudinal deformations. Comparing 
the simulated MR image and the real follow-up of the success¬ 
ful case, the tendency for the ventricles to increase in size, 
reflecting more cortical atrophy was successfully simulated. 
The VBM map also showed an approximately matched tissue 
deformation along the ventricle as well as the cortical areas. 
The example that failed showed a mismatch of the atrophy 
in the MR images and the VBM maps. The increase in 
size of the lateral ventricles was neglected and the cortical 
atrophy was overestimated. Such simulation failure can be 
generally attributed to the affine normalisation errors and the 
insufficiency of the longitudinal template collection. 

To evaluate the strengths of different template weighting 
strategies, we considered the square distances between the 
simulated volume and the real follow-up volume (P-B) as 
well as the distances between the simulated volume and the 
registered Year 1 volume (P-rA), which is the output of 
registering the Year 1 image to its real follow-up MR image 
(Fig. [2]). In Fig. [2] unlike P-rA, the P-B curve is correlated with 
the original registration error. It might indicate that P-rA could 
be a more unbiased evaluation criterion for such simulations. 
All the distances in Fig. [3] are sorted according to the image 
distances of the voxel intensity based simulation (intensity). 
The distances were zero-meaned and rescaled to make the 
individual differences identifiable. In Fig. [3j we compared the 


longitudinal (long) and cross-sectional (cross) simulations as 
well as the combined weights of both (combine) respectively 
according to the P-B distances (Fig. |3]-(a)) and the P-rA 
distances (Fig. |3]-(b)). It is noticeable that the longitudinal in¬ 
formation (long) outperformed the cross-sectional information 
(cross) in most cases. Since neither method achieved the lowest 
simulation errors in all trials, the combined weights (combine) 
could be used to balance two perspectives. Taking the intensity- 
based method (intensity) as a reference, at least one of the 
proposed morphometry-based methods (long, cross, combine) 
achieved lower simulation errors in most cases regarding both 
criteria. 


IV. Conclusion 

We present a framework which automatically simulates fu¬ 
ture neurodegeneration from the longitudinal atrophic changes 
as well as the cross-sectional difference from a statistical 
average normal template. This framework expects the patient to 
have at least 2 historical MR images to increase the confidence 
of the simulated 3D MR volume. The brain tissue deformation 
was represented by the voxel-based morphometry (VBM). Our 
evaluation showed that the intra-subject longitudinal infor¬ 
mation enhances the simulation accuracy. The results from 
at least one of the proposed morphometry-based methods 
outperformed the state-of-the-art intensity-based method in 
most evaluated cases. With a sufficient template collection, our 
proposed framework can be used for validating the prediction 
made by neuroimaging measurements extracted from MR data. 
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(a) Sorted distances between the predicted MR and the real follow-up MR (P-B) 



(b) Sorted distances between the predicted and the registered Year 1 (P-rA) 


Fig. 3. The plot of different template weighting methods: cross-sectional VBM weighting (cross), longitudinal VBM weighting (long), combined weighting 
(combine), intensity-based weighting (intensity). All distances were zero-mean rescaled into [—1,1] and the subjects were sorted according to the intensity-based 
weighting. Lower values indicate better simulation results. 
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